clear all; close all; clc;

H = 19100e3;
Rz = 6400e3;
Tobr = 12 *60 *60;
Tmod = Tobr /2;
wobr = 2 *pi / Tobr;
c_light = 3e8;
f_0 = 1602e6;

dt = 5*60;
t = 0:dt:Tmod;

rr=[0; Rz];

D = nan(1, length(t));

for i=1:length(t)
   rsv = (H + Rz)*[cos(wobr*t(i));sin(wobr*t(i))];
%    plot([0 rsv(1)],[0 rsv(2)]);
%    xlim([-H-Rz; H+Rz])
%    ylim([-Rz; H+Rz])
%    axis equal
%    pause(0.1);
   rv = rsv -rr;
   
   D(i) = norm(rv);% дальность
   y_rv(i)=rv(2);
end

%получаем скорость из дальности
%функция diff сделает массив на 1 элемент меньше

ind = find(y_rv > 0);
V = diff(D)/dt;

df= - V/ c_light* f_0;
figure; 
plot(t(2:end),df, t(ind), df(ind-1));
xlabel('t, sec')
ylabel('\Delta f, Hz')

figure
plot(y_rv)

figure;
hist(df)
xlabel('Delta f, Hz')